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Abstract This paper presents the general spatial Schmidt decomposition of two-photon fields generated 
in spontaneous parametric down- conversion (SPDC). It discusses in particular the separation of the radial 
and azimuthal degrees of freedom, the role of projection in modal analysis, and the benefits of collinear 
phase mismatch. The paper is written in a review style and presents a wealth of numerical results. It aims 
at emphasising the physics beyond the mathematics, through discussions and graphical representations of 
key results. The two main conclusions of the paper are the finding of a better law to describe the effective 
dimensionality of the spatial part of the total Hilbert space and a possible novel feature of the radial 
Schmidt modes. 
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1 Introduction 

The understanding of processes that generate entangled 
states is of great importance. Current research makes a 
vast use of entangled states, and future technologies might 
rely on quantum properties of matter and fields [1 . Spon- 
taneous Parametric Down-conversion (SPDC) is one of 
such processes, which can generate high-dimensionally en- 
tangled states of two or more photons [2 l3|4|5|6j . The lit- 
erature on SPDC is well developed, however there are still 
some open questions, some of which we aim at answer- 
ing. Although the entanglement between SPDC photons 
happens also between degrees of freedom beyond the spa- 
tial ones (polarisation and time- frequency for instance [7 ), 
given the type of experimental setups that we refer to, we 
will concentrate on the spatial part of the entanglement. 
Therefore, the key challenge we will address in this paper 
is how to decompose the two-photon field A{cii^q^2) into 
its Schmidt modes lii(qi) and Vi{q^2)- 



^(qi,q2) 



E 



AiUi(qi)wi(q2), 



(1) 



In particular, the issue of being able to take into account 
a phase mismatch between the photons, when calculat- 
ing entanglement strength and detection probabilities has 
not yet found a full resolution. Approximations have been 
made, and with these, very interesting analytical results 
have been found [6 8 . However, approximations naturally 
restrict the domain of validity of the results. The most 
common approximation used in SPDC theory (replacing 
a sine function with a Gaussian function) is discussed in 
this paper. We present a numerical analysis of the en- 
tanglement strength, of the Schmidt modes, and of the 
detection amplitudes without the usual approximations. 
The first part of the paper is a theoretical discussion of 



the physical system under study and a presentation of the 
meaning of the Schmidt decomposition. In the second part 
we present the numerical results, and we describe some in- 
teresting features of the structure of the entanglement. 



2 Theory 

2.1 Generated two-photon field 

We consider the generation of entangled photon pairs via 
the nonlinear optical process of SPDC, where single pump 
photons which propagate through a nonlinear crystal oc- 
casionally split in two photons of lower energy, tradition- 
ally called the signal and idler photons. We operate in 
the quasi-monochromatic frequency-degenerate limit, by 
combining a narrow-band cw pump laser at frequency ujp 
with narrow-band spectral filters at the degenerate fre- 
quency ui = 002 = ujp/2. When the pump field is spa- 
tially coherent and the spectral filtering is sufficiently nar- 
row, the generated two-photon state (ignoring polarisation 
and therefore describing both type-I and type-II down- 
conversion) is pure and has the general form [9] 



(2) 



where d^qj^) is the creation operator of a plane wave with 
transverse momentum and i is the photon label. 

One of the factors of the two-photon amplitude A(qi, q2) 
in eq. ([2|, is the phase matching function, which describes 
the efficiency of the down- conversion process, as a function 
of transverse k-vectors. The longitudinal phase matching 
function in the two-photon amplitude is a sine function, 
being the Fourier transform of the step function that rep- 
resents the uniform amplitude of creating the photon pair 
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along the crystal. It is usually approximated with a Gaus- 
sian function if no or very little phase mismatch (indi- 
cated in the paper by (p) is present. However, the descrip- 
tion that derives from such approximation cannot account 
for an arbitrary phase mismatch. Moreover, even in the 
regime where ^ = the difference in shape of the approx- 
imated amplitude from its complete form has an effect in 
the analysis of the down converted state, as we will show. 

The the two-photon amplitude in the far field, or in 
momentum space, presents a ringed pattern, due to the 
sine function. If the phase is matched ((/? = 0) the sine 
gives a bright spot in the centre and secondary rings around 
it. Most of the intensity of the field is concentrated in the 
central area, it is for this reason that in case of perfect 
phase matching it is possible to roughly approximate the 
phase matching function with a gaussian function (see be- 
low for more details). In this paper we also analyse situa- 
tions in which ^ 0. In this case the argument of the sine 
function is modified by a constant value, see eq. 0. The 
effect in the amplitude for negative and decreasing ip is 
to lose the central spot and obtain a more divergent field, 
in fact every secondary ring obtains a larger radius the 
lower the value of Lp. On the other hand, for positive and 
increasing ip the effect is to obtain an overall weaker field, 
less divergent, with a central spot intensity that depends 
on the value of Lp. We will refer to the effect of tuning Lp 
as opening and closing the rings. 



2.1.1 Sine phase matching 

The complete form of the two photon amplitude is [10 

A(qi, q2) (X fp(qi + ^2) X sine [\AkzL) , (3) 

where fp(qi + q2) is the transverse momentum profile of 
the pump beam (q^ = qi+q2). The function sine {^AkzL) 
quantifies the influence of phase matching, where L is the 
crystal length and where the projected wave-vector mis- 
match Akz = Akz{c[i,q^2) is a function of qi and q2. 

We make the general Eq. (|3| more specific by inserting 
a (rotationally-symmetric) Gaussian pump profile. Other, 
more general, types of pump profiles can be treated analyt- 
ically in the limit for a short crystal [11 . We also assume 
non-critical phase matching, to remove any linear depen- 
dence of Akz on Qi? and perform a Taylor expansion of 
^^2(^1,^2), neglecting a small term that scales with the 
difference between the refractive indices at the pump and 
SPDC wavelength multiplied by |qi + q2p- Following the 
notation of Law and Eberly [10 , we thus rewrite Eq. ([3| 

A(qi,q2) (xexp(-|qi + q2 1 V^^^) xsinc (6^|qi -q2p + V^) , 

(4) 

where a = for a Gaussian pump profile fp(x) cx 

exp (— |xp/i(;^), and 6^ = L/{Akp)^ where kp = nup/c is 
the momentum of a pump photon in the crystal [12 , and 
(p is the collinear phase mismatch. 

Although this challenge has been tackled mathemati- 
cally by Law and Eberly [10], several aspects and phys- 
ical implications of the Schmidt decomposition remained 



undiscussed and are still a topic of active research |6|8j . 
This paper aims to be as complete as possible in this dis- 
cussion. 



2.1.2 Gaussian phase matching 

As anticipated, to obtain analytic solutions, many authors 
|6|10|12j replace the sinc-type phase matching function at 
(p = by its Gaussian approximation, using the substitu- 
tion sine (6^|qi — q2p) ^ exp (— 6^|qi - q2p). The Gaus- 
sian approximation allows an exact Schmidt decomposi- 
tion with Hermite- Gaussian or Laguerre- Gaussian eigen- 
modes |6|10|8|13] . The Schmidt number of this decompo- 
sition, i.e., the eflFective number of modes that participate 
in the modal decomposition, is [10 

K=\[^^l)\ ,5) 

Most experiments operate in the weak- focusing limit (large 
Wp^ long Rayleigh range), where 6cr <C 1 and K ^ l/{2ba)'^. 
The generic beam waist of all Hermite-Gauss (HG) or 
Laguerre- Gauss (LG) modes, defined by the fundamen- 
tal mode iXo,o(q) = exp (— g^/g'Q), is qq = y^a/{2b) being 
the geometric mean between the width of the pump profile 
and the Gaussian phase matching function. 

The eigenmodes of the Gaussian Schmidt decomposi- 
tion are not unique (because of the degeneracy of all the 
eigenmodes belonging to a fixed mode order) and, for ex- 
ample, can equally well be taken as the Hermite- Gaussian 
or Laguerre-Gaussian set of modes. 

The numerical results presented in section 3 suggest 
that the Gaussian expansion can be improved somewhat 
by writing sine (6^|qi - q2p) ~ exp (-a^6^|qi - q2p), 
where a is a (scaling) constant chosen such that both func- 
tions satisfy some common criterium. This modified scal- 
ing results in a replacement of b ^ ab in all expressions 
that originate from the Gaussian approximation. We post- 
pone the discussion on the rescaling factor a at the end of 
the paper, as additional concepts have to be introduced. 
In Figure [1] we plot the sine phase matching (solid green 
line), comparing it to its Gaussian (red dashed) and super- 
gaussian (blue dashed) approximations. It is straightfor- 
ward to notice that the Gaussian approximation exp(— g'^) 
scales very rapidly, thus missing out high-g modes. The 
sine function sinc(g'^), instead, has di overall scaling, 
letting high-g' modes have an effect on the Schmidt decom- 
position. The supergaussian exp{—bq'^) with b = 0.193 is 
a better approximation than the Gaussian, as it captures 
the central scaling of the sine much better |14|15j . This 
is due to the Taylor expansion of the supergaussian hav- 
ing the same powers as the sine. However, there are two 
main problems associated with the supergaussian approx- 
imation. The first is that although it models the central 
scaling for a wider range of q values, it misses the 1/q^ 
overall scaling too, failing to carry high-g' modes, and it 
does not display any ring structure. The second problem is 
that analytically it is not an advantage: for perfect phase 
matching it is possible to find angular Schmidt modes in 



Filippo M. Miatto et al.: Spatial Schmidt modes generated in parametric down-conversion 



3 



terms of Bessel functions of two variables, which are de- 
fined in [16 17 , but as the literature on generalised Bessel 
functions is poorer than the literature on single variable 
Bessel functions, it is perhaps an unprofitable effort to 
try to find the full analytical Schmidt decomposition us- 
ing the supergaussian approximation, and as we will show 
in section 3.3, a heuristic rescaling of the Gaussian ap- 
proximation is sufficient to express the Schmidt number 
correctly in case of perfect phase matching. 




Figure 1. (Color online) Comparison between the three 
different phase matching functions. The solid line is the 
full sine form sinc(g'^), the dashed lines are the gaus- 
sian (red, narrower) exp(— g^) and optimised supergaus- 
sian (blue) approximation exp(— 0.193g'^). The similarity 
between the Gaussian and sine functions is questionable, 
as it is the similarity between supergaussan and sine for 
the overall scaling and for the oscillations at <7 ^ 2. 



2.2 Schmidt decomposition and Entanglement 

For the arguments in this paper it is sufficient to recall that 
the Schmidt decomposition of an entangled state is its rep- 
resentation in the joint basis that maximises the mutual 
information. Notice that in Eq. ([T]) the joint modes are the 
product of one mode per subsystem (hence each Schmidt 
mode is separable). The decomposition, in fact, also speci- 
fies the weights of single joint modes '^i(qi)'^i(q2), the 
square of which is interpreted as the normalised probabil- 
ity (via = 1) of detecting the i-th joint mode, were 
we to measure in the Schmidt basis. The more the distri- 
bution of the probabilities {Xi} is "spread out" the more 
entangled the state is. This leads to the natural choice of 
choosing entropic measures to quantify entanglement. In 
our case we choose the Schmidt number K = l/X^i^i^ 
where the probabilities are normalised: = 1. 

The Schmidt decomposition of the state ([2| can be 
thought as made of two stages, because we are working 
with two pairs of coordinates: an angular coordinate per 
subsystem and a radial coordinate per subsystem. In the 
first stage we separate the angular variables, in the second 
stage we separate the radial variables. The decomposition 



will have therefore two sums, and the weights and joint 
modes will have two indices [6 . 

We start with the two-photon amplitude 74(qi, q2) and 
we write it in polar coordinates (which we call qi^2 and 
0i^2)' The rotational symmetry of A(qi, q2) limits the an- 
gular dependance to a function of {6i — O2) only. The or- 
thogonality of the azimuthal functions e*^^^^"^^) allows for 
an easy extraction of the quantity 




(91,92)= / A{ciuCi2)e-'"^''-'^U{9,-92). (6) 



^0 

This means that with a Fourier transform we can write 
the two photon amplitude in a joint basis of the Orbital 
Angular Momentum (0AM): 

^(qi,q2) = ^ E VPtmi,<i2)e''''e-^''^ (7) 

Where is the probability of measuring a pair of pho- 
tons with 0AM ^ and — ^, were we to measure in the 
joint 0AM basis. The radial amplitude is normalised as 
J \Fi{Qi^Q2)\'^QiQ2dQidQ2 = 1- In this way we have sep- 
arated the two angular variables, and each joint 0AM 
mode is labelled by an integer number £. We can therefore 
proceed to the second stage and consider the functions 
\^F£{qi^q2) and separate the radial variables. In order 
to make this separation easy and truly one-dimensional, 
we follow the reasoning of Law and Eberly [10^ by intro- 
ducing a prefactor y^^i^, to account for the Jacobian of 
the transformation q and separating the radial 

function as 

^/q^\fP^tFi{qi,q2) = ^ \/A^(/>£,p(gi)0-£,p(g2). (8) 

Symmetry considerations also lead to the separated func- 
tions to be equal, i.e. we use (j){qi) as well as (/>(<72)- The 
prefactor yjq\q2 is needed to facilitate the simple nor- 
malisation of the eigenfunctions (J (l)i^p{q)(l)i^p'{q)dq = 
Sp^p') that is required for the numerical decomposition that 
we discuss later on. For convenience, we will also define 
^Ap(<?) = (t)£,p{q)/^/q such that f Ui^p{q)ui^p>{q)qdq = 8p^p>. 

Note that we had to introduce the non negative, in- 
teger label p. Such quantum number labels radial eigen- 
modes. The numerical analysis will lead to some more in- 
sight on this quantum number, as we will show in sec- 
tion [3] After putting these two results together we obtain 
the full two-dimensional Schmidt form of the two-photon 
amplitude 

00 CO 

^(qi,q2)= ^ \/A^^Ap(qi)^-Ap(q2). (9) 

Note that this form is analogous to Eq. ([T]) with Ui{c\) = 

ui^p{q)e'^^. 

The radial shape (l)£^p{q) of these profiles depends stron- 
gly on the collinear phase mismatch cp; it resembles the 
well-known Laguerre- Gaussian (^,p)-profiles only for (f ^ 
0, but deviates strongly from these profiles at ^ (see 
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below). As anticipated, the effective number of spatial 
modes (which can be thought of as the strength of the 
entanglement) is given by the Schmidt number 



K 



(10) 



where the probabilities are normalised: Xl^p ~ ^^ 

A simple scaling argument shows that the outcome of 
the decomposition depends only on two dimensionless pa- 
rameters. The first parameter is the product ba = . / ^^^^1^ = 



^/L^J2 where we define to be the crystal thickness 
normalised to the Rayleigh range of the pump beam, = 
L/zr. The second parameter is the collinear phase mis- 
match ip. 



2.2.1 Entanglement in orbital angular momentum 

The rotational symmetry of the considered geometry and 
the associated conservation of 0AM makes the azimuthal 
component of the entanglement easier to address and more 
fundamental than its radial counterpart. Hence, most ex- 
periments have concentrated on the 0AM part of the spa- 
tial entanglement. In this context we are interested in the 
probability Pi (which appears in Eq. ^) of detecting a 
pair of photons with 0AM, respectively, £ and —£. Such 
probability is related to the weights Xi^p via Pi = A^^p. 
Pi is called the spiral weighty its full distribution {P^} is 
called the spiral spectrum^ and the width of this distribu- 
tion is called the spiral bandwidth [3 . The dimensionality 
of the generated 0AM entanglement is defined by the az- 
imuthal Schmidt number 



(11) 



To be compared to the full Schmidt number ([To]). In case 
of high-dimensional spatial entanglement {K ^ 1) and 
(p ^ 0, the relation between both forms of entanglement is 
i^az ^ 2V^, as the effective range of the £ labels is typi- 
cally 4x as large as that of the p labels. This statement 
is exact for the Gaussian approximation of the two-photon 
field, where [13] 



Xe,p = C,exp{-{2p+\£\)/C2), 



(12) 



where Ci and C2 are two constants. It does, however, not 
apply to the full expression for cp ^ 0. Our numerical 
approach allows us to analyse what happens also in such 
regime. 



2.3 Optical etendue 

The Schmidt decomposition of the two-photon field is equiv- 
alent to a coherent mode decomposition of the coherence 
function of the reduced one-photon state. The reduced 
one-photon operator obtained by tracing |^)(^| over 
all possible states of the second photon, is an incoherent 



(and weighted) mixture of all of all Schmidt states. The 
extent to which the reduced density operator is mixed re- 
flects the degree of entanglement for the pure two-photon 
state IV^). The Schmidt number K can thus be estimated 
by comparing the one-photon coherence of the SPDC emis- 
sion with that of more standard optical sources. 

At sufficiently large ^ 1, the reduced one-photon 
field is quasi-homogeneous, i.e., the reduced one-photon 
density matrix p*^^^(qi,qi) factorizes in a product of an 
intensity function of sum coordinate (qi -f q'^^) and a co- 
herence function of difference coordinate (qi — q'^). The 
Schmidt number K of the quantum state then reduces to 
the normalised optical etendue N = Af2/X^ of the one- 
photon field, where A is the effective area of the source 
and i? is its effective (solid) opening angle [13 18 . The 
optical etendue quantifies the effective number of trans- 
verse modes contained in a partially-coherent beam. The 
relation K ^ N becomes exact if the source is quasi- 
homogeneous and if the effective area A and opening angle 
i? are defined in a convenient way [18]. 

The relation between K and N yields an easy and in- 
tuitive interpretation of the Schmidt dimension. For this, 
we convert the opening angle of the source i? into a trans- 
verse coherence width Wcoh in the source plane, which we 
conveniently define as 



^coh = 46 = 



nn 



(13) 



where L is the crystal length, Aq is the emission wave- 
length in vacuum and n is the refractive index. The trans- 
verse coherence width Wcoh of the one-photon field is sim- 
ilar to the width of the two-photon coherence function 
V{xi—X2) defined in refs. jl9|1 8 as the Fourier-transform 
of the sine- type phase matching function. With the above 
definition of Wcoh^ the Law-and-Eberly expression for the 
Schmidt number at 6cr <C 1 reduces to the logical form 



K 



1 

2b^ 



Wcoh J 



(14) 



which shows that the Schmidt number K simply counts 
the number of "independent" coherent regions in the source. 
The above definition of Wcoh also provides for an easy 
rewrite of the results of the Gaussian expansion discussed 
in subsection (2.1.2). The width wq of the fundamental 
Gaussian Schmidt mode uo^o{q) oc exp— in real 
space is Wq = 2/qo = y/WpW coh^ being the geometric 
mean between the pump waist and the transverse coher- 
ence length of the generated field. The modal amplitude 
of the Gaussian U£^p mode is 



exp[-{2p+\i\)/KiD] 



(15) 



where Kuj = Wp/wcoh is the one-dimensional equivalent 
of the Schmidt number. 



2.4 Projective measurement of 0AM entanglement 

In this subsection, we will briefly compare two experimen- 
tal techniques that have been developed to characterise in 
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particular the 0AM contents of two-photon sources. We 
distinguish between measurements with bucket detectors, 
which record the complete field, and measurements with 
single-mode detectors, which record projected components 
of the field. 

A Schmidt analysis of the two-photon field with bucket 
detectors, which by definition have no spatial sensitiv- 
ity, requires two-photon interference before measurement. 
More specifically, one should measure the visibility of the 
famous two-photon (Hong-Ou-Mandel) dip [20 as a func- 
tion of the relative alignment of the interfering beams. 
Measurements with one beam rotated with respect to the 
other provide complete and unbiased information on the 
azimuthal part of the entanglement |21|22j , as the visibil- 
ity of the observed two-photon interference is 

V{Ae)^ J J A*{qi,q2;A0)A{qi,q2;A0 + A0)qiq2dqidq2 



Peexp{iiAe), 



(16) 



where is the rotation angle between the interfering two- 
photon fields, which is twice the rotation angle applied 
in a single arm of the interferometer. This angular depen- 
dence originates from the interference between two-photon 
fields with either the i or the —£ photon in arm 1 and the 
other photon in arm 2. Equation (16) is based on the en- 



tanglement between these contributions and the natural 
symmetry P_£ = P^. 

A Schmidt analysis of the two-photon field with single- 
mode detectors, each comprising a single-mode fiber posi- 
tioned in front of a detector, requires some form of mode 
transformation before projection. This mode transforma- 
tion is typically performed with (a set of) fixed phase 
plates or an (adjustable) spatial light modulator (SLM). 
The combined transformation-projection technique has been 
used successfully for the analysis of the 0AM contents of 
the two-photon field [18 . The generated two-photon field 
is then typically projected onto two modes of the form 
(l)d{q)ex.-p{if{6)), where 0d(g) is the (image of the) mode 
profile of the fiber and exp {if{0)) is the phase profile im- 
posed by the phase-transforming element, which typically 
has the standard form exp {i£0) . The rotational symmetry 
of the generated field 74(qi,q2), which manifests itself in 
the generation of photon pairs with opposite 0AM only, 
allows one to derive simple expressions for the projected 
two-photon field. It is easy to show that the maximum 
information contained in these projections, for any com- 
bination of azimuthal phase profiles /(^), is contained in 
the function 

^proj(/A(9) = J J M^i^^2)K{qi)u'^{q2)qiq2dqidq2 

oo oo 

cx ^ ^Q,pV^exp(iM(9), (17) 

where AO = Oi — O2. The projection coefficients 

Ci^p = I / Ui^p{q)ud{q)qdq\'^ < 1, (18) 



quantify the spatial overlap between each generated Schmidt 
mode and detector mode. They also quantify the bias im- 
posed by the mode projections. 

As an example, let us consider the projection on a 
detector mode Ud{q) that is much wider than a range 
of lower-order U£^p{q) modes. Under these conditions, the 
projection coefficients will have a strong bias towards the 
p = modes, as all U£^p{q) with p > 1 modes exhibit 
at least one oscillation (in the radial direction) making 
Q,p/o <^ Q,p=o- Hence, projection experiments are mainly 
sensitive to the p = components of the two-photon 
field. A calculation of the projection-imposed bias on the i 
modes is more complicated. A rigorous analysis, in terms 
of Lerch transcendent functions, is presented in another 
contribution to this issue 23l. 



3 Numerical results 

In this section we will presents numerical results for the 
general Schmidt decomposition of Eq. Q into Eq. ([9|. 
We will in particular discuss the distribution of the modal 
weights Xi^p, and their separation over the azimuthal and 
radial degrees of freedom, the shape of the Schmidt modes 
U£^p{q)^ and the influence of the collinear phase mismatch 
(p. 

The numerical analysis is relatively straightforward. 
We start from the generated two-photon field 74(qi,q2) 
defined by Eq. Q and perform an 0AM decomposition 
into radial-only function V^^^(^i7^2) using Eq. The 
analytical integration is performed using an approxima- 
tion that works well for bcr < 0.2, which is the most rele- 
vant experimental situation as it corresponds to a weakly 
focused pump. Employing the analytical integration of the 
angular variable yields analogous final results in much less 
time, as the numerical methods consist only in a diagonal- 
isation of matrices. The prescribed radial integration is 
performed on a discrete equidistant grid of g'- values. After 
multiplication by a factor ^Jq\q2^ to account for the Jaco- 
bian of the transformation, we obtain the scaled N x N 
matrix F with a dimension set by the required accuracy 
{N is typically in the order of the hundreds) . The Schmidt 
decomposition in Eq. ^ of \fPiFi{qi^q2) is equivalent to 
the diagonalisation of the matrix F. To avoid potential 
complications associated with the two-particle nature of 
the ^(^1)^(^2) phase space, we first multiply F by its 
transpose F''', in order to obtain the equivalent of the 
reduced one-photon density matrix, and diagonalize the 
matrix FF^^ = F^ instead. The diagonalisation of the ma- 
trix yields three matrices: = MEM^ ^ two of which 
represent the modes (M) and one of which is diagonal 
{E) and represents the eigenvalues. The resulting eigen- 
values are X^^p (at fixed ^); the resulting eigenmodes are 
</>Ap = V^Ui^p{q). 

In the following subsection we present the weights in 
the matrix E and we show what are the consequences 
of having a nonzero phase mismatch. Similarly, we also 
present the modes in the matrix M, and we show the 
consequences of having a nonzero phase mismatch. 
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3.1 Schmidt weights 

As a typical example we will concentrate on the case ba = 
0.05, which for instance corresponds to a 2 mm thick PP- 
KTP pumped hy di Wp = 230 fim diameter pump laser at 
a pump wavelength of 413 nm. 



3.1.1 Phase mismatch (p = 
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Figure 2. False color plot of the modal spectrum Xi^p 
(with p horizontal and £ vertical) for ba = 0.05 and (p = 0. 
Note the typical triangular form and the ^ 4x larger range 
of relevant £ values as compared to p- values. 



Figure [2] is a false-color plot of the probabilities X^^p of 
the modal aecomposit ion for ba = 0.05 and (p = 0. The 
Schmidt numbers deduced from this calculation are K ^ 
231 and iCaz ^ 32, respectively. These numbers satisfy 
the approximate relation i^az ^ 2^/K associated with the 
triangular mode spectrum in figure [2] They are, however, 
considerably larger than the values of K ^ 100 and Kg,z ^ 
20 expected from a simple Gaussian expansion based on 
the approximation sinc(g^) ^ exp(— g'^). The reason for 
this difference, as anticipated, is the fact that the Gaussian 
approximation is missing high-g transverse modes. 



3.1.2 Phase mismatch (p < 

We consider now geometries with non-perfect phase match- 
ing and in particular the case (p < for which the SPDC 
rings have opened up. 

Figure |3] is a false colour plot of the probability X^^p of 
the modal oecomposition calculated for ba = 0.05 and (p = 
—4. On account of phase mismatch, the full 2-dimensional 
Schmidt number has increased by a factor ^ 2x to K ^ 
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Figure 3. False color plot of modal spectrum A^^p (with 
p horizontal and £ vertical) for ba = 0.05 and cp = —4. 
Note the strongly elongated form in the £ direction and 
the reduced width in p. 



425. This increase coincides with a similar increase in 
the space angle and the angular integrated output of the 
SPDC source. Most noticeably, the spread in 0AM values 
is now much larger then in the (p = case, at the expense 
of the spread in p values. The azimuthal Schmidt num- 
ber has increased from its prior value of 32 at = to 

i^az ^ 72. 



3.1.3 Phase mismatch (p > 

For the sake of completeness we consider now the case of 
positive phase mismatch, (/? > 0. In this case the SPDC 
rings have closed, and the efficiency of the down-conversion 
process is considerably lower, normally by two orders of 
magnitude, which degrades the S/N ratio by a large ex- 
tent. 

Figure [4] is a false color plot of the eigenvalue Xi^p of 
the modal aecomposition calculated for ba — 0.05 and 
= 4. As compared to the case of negative detuning, we 
find that the width in the p direction has increased at the 
expense of the width in the ^-direction. In this example 
we calculated values oi K 625 and Kaz ^57. 



The optical etendue argument presented in Sec. (2.3) 



allows one to explain why the spatial entanglement gets 
concentrated in the azimuthal degree of freedom, instead 
of the radial one, iox ip < and why this concentra- 
tion hardly occurs for > 0. Our geometric argument is 
based on the notion that the transverse coherence width in 
the azimuthal and radial directions are comparable, being 
Fourier related to the size of the pump beam. At negative 
detuning (p <0^ where the fundamental SPDC ring opens 
up, the circumference of this SPDC increases with phase 
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Figure 5. Probability distributions of 0AM modes (i.e. 
summing over p) . The Lorenzian-like distribution (red) is 
for positive mismatch and the Gaussian-like distribution 
(blue) is for negative phase mismatch. 



Figure 4. False color plot of modal spectrum A^^p (with 
p horizontal and i vertical) for ba = 0.05 and cp = 4. 
Note the strongly elongated form in the p direction and 
the reduced width in £ as compared to the case cp = —4. 



mismatch while the radial thickness of the ring decrease 
accordingly. As a result, the number of spatial modes that 
fit in the azimuthal {£) direction, increase while the num- 
ber of modes in the radial {p) direction decrease. The to- 
tal area of the open SPDC ring at <C — 1 is approxi- 
mately twice as large as the area of central emission at 
(/? = 0, making the 2-dimensional K also about two times 
larger. At positive detuning (p > the fundamental SPDC 
ring disappears and only weak secondary rings remain. As 
these rings are numerous and have similar intensities, the 
number of spatial modes in the Schmidt decomposition 
can be quite large and is more evenly distributed over £ 
and p mode numbers. 

Next we single out the 0AM part of the entanglement 
for the considered cases (/p = 0, — 4 and 4. We do so by sum- 
ming the modal weights A^^p, depicted in figures |2j [3] and [4| 
over its radial quantum number p. What one then obtains 
are three different probability distributions that show dif- 
ferent behaviour over the same range of £ values. Figure 
Is] shows the modal weights (^^ ^i,p) of the 0AM modes. 
Note that the behaviour changes from Lorenzian-like to 
Gaussian-like as the phase mismatch goes from positive 
to negative values p2]. We stress that the case of posi- 
tive detuning is experimentally unfavorable on account of 
the limited brightness of a source where the fundamental 
SPDC ring has closed/disappeared. 



3.2 Schmidt modes 

In this subsection we will present a series of numerical 
results that show the radial profile of the Schmidt modes 
in momentum space. The profile of the modes is influenced 



by the presence of a phase mismatch, and as anticipated 
the modes resemble the LG modes only at = 0. 



3.2.1 Phase mismatch cp = 

Figure |6] shows the result of the singular value decompo- 
sition for ba — 0.05. What is represented is a portion of 
the matrix M corresponding to a value of the 0AM of 10. 
This choice gives a good example of the radial displace- 
ment in /c-space between Schmidt modes that we observe. 
Such behaviour will become clear below, where we com- 
pare modes from a perfect phase matched case to modes 
from a negatively detuned phase mismatch. Each horizon- 
tal line represents the intensity profile |(/)£^p((7)p of a radial 
mode 0^,p(^) with fixed £ and a different radial quantum 
number p. The first three modes in figure |6] are plotted in 
figure [71 

We would like to give an interpretation of the be- 
haviour of the Schmidt modes in figure [6] Notice that 
each eigenmode resembles a common eigenmode from the 
family of Hermite-Gauss modes, only centred at specific 
distances from the origin. There is a one to one correspon- 
dence between a family of modes and a region between ze- 
ros in the function \/^^^(^i7 ^2)- Such function represents 
the joint amplitude density of the photons having radial 
momentum qx _and ^2, respectively (examples are given in 
figures 



12 



and 13). Therefore it's also closely related to the 
radial sine profile of the down- converted beam. Such cor- 
respondence might be explained by considering the size of 
the coherence length compared to the size of the nonzero 
regions in V^^^(<?i7 ^2)- if the coherence length is smaller 
than the thickness of a ring, there can be no strong co- 
herence between photons in different rings, and any wave 
function in a single ring can be written in terms of a set 
of fundamental modes of fixed 0AM and variable p. 

We find the organisation of the eigenmodes in groups 
quite remarkable. Each of the groups looks like a local 
complete family of HG-like modes, giving the impression 
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of the existence of an additional symmetry and therefore 
of a quantum number that could be specified to identify 
each group. This can be explained if the matrix is taken 
into account and in particular the shape of the lobes that 



0.25 



form along the diagonal, as in figure 12 and 13 They can 
be thought of as potential wells, which have to be filled 
by eigenmodes. For each of them the eigenmodes build in 
local families because from the Taylor expansion of the 
neighbourhood of the bottom of any potential well we can 
infer that the HG modes form a natural set of local modes 
for each well. The potential walls between neighbouring 
wells could be tunneled by the modes, if their size were 
large enough to reach over to the next well. In practice, 
the size of each mode is small enough and the families are 
well confined within individual lobes of the sinc-fucntion 
if the pump beam is not too focussed. 
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k- space radial variable q 

Figure 6. A part of the matrix M for ^ = 10 is plotted in a 
black and white array plot. The intensity profile 
of a single radial eigenmode is represented in each hor- 
izontal series. The radial quantum number p labels the 
different series. They are organised from top to bottom in 
order of contribution to the total wave function. In this 
example the radial field that was decomposed was the one 
in figure [12] 
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Figure 7. Amplitudes in the form 
three p- modes in figure [6] 



]^p{q) of the first 



As with the zero mismatch example, we also plot the first 
three modes of figure [8] in figure |9] Note that although 
the nonzero phase mismatch, the shape of the first three 
eigenmodes is remarkably similar to the one in figure [7| the 
difference being just an offset with respect to the origin. 
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Figure 8. A part of the matrix M for ^ = in case of 
phase mismatch (p = —4. In this example the radial field 
that was decomposed was the one in figure [13) 



3.2.2 Phase mismatch (p 

Figure [s] shows the calculated radial profile \(j)£^p{q)\'^ of 
the eigenmodes at £ = 10 for the detuned case (p = —4. 
At (/9 <C 0, where the SPDC ring has opened, even the 
fundamental i = p = Schmidt mode is ring-shaped and 
has a wide region of zero intensity around the central axis. 



This happens because in the two cases the radial func- 
tions that are being decomposed are remarkably similar. 
The modes firstly fill the largest ring, until the difference 
between their combined intensity and the intensity of the 
main ring is small enough to be less than the intensity 
of secondary rings. When this happens we see that the 
subsequent modes fill different rings, and jump discretely 
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between them. This is in accordance with the explanation 
that we gave at zero phase mismatch: the modes are com- 
pact enough to be well contained within individual rings. 



0.25 
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Figure 9. Amplitudes in the form (j)i^p{q) of the first 
three p- modes in figure |8] 



3.3 Heuristic rescaling of Gaussian approximation 

In this section we conclude the argument, anticipated in 
the previous sections, on the scaling factor a. Such scal- 
ing factor would allow us to describe more accurately, for 
(/9 = 0, the Schmidt number when approximating the sine 
phase matching with a Gaussian phase matching. Below, 
we present two arguments for why we think it is better 
to replace the function s\nc(}p'q^) by exp {—a^h'^q^) rather 
than by exp {—Ij^q^). There are experimental results that 
support our claims, for instance see [24 . 

Figure [lO] depicts the numerically-calculated Schmidt 
number K and its azimuthal counterpart i^Taz as func- 
tions of ha for the exact sine phase matching. In order 
to numerically calculate these curves we had to abandon 
momentarily the faster approach used to find the Schmidt 



modes, as her in figure 10 is not restricted to small values. 
The Gaussian approximation of the phase matching func- 
tion results in expressions that are symmetric under the 
exchange ha ^ The same argument cannot be applied 
to the sine phase matching, as it's possible to see from 
figure [To] the graph is symmetric with respect to a point 
between 1 and 2, meaning that there could be an extra 
factor a such that the symmetry now reads haa ^ 
representing a reflection around the point 1/a. If it were 
so, the factor a should have a value around a o:^ 0.85. 
However, the introduction of a alone is not enough to ob- 
tain the correct K values. There are two ways to proceed, 
either to choose a simple modification, like ha haa^ or 
to add also a second multiplicative factor p. We find that 
the relation 



K=^ (-^^haa 
4 \haa 



fits extremely well the values in figure [To] for a = 0.85 
and P = 1.65. However, if one is only interested in the 
regime 6cr < 1, which experimentally is the most sig- 
nicative, there is no need to introduce /3 and the relation 
K ^ l/{2haa')'^ works well for a' = a/y^ ^ 0.65, as it is 
possible to see in figure 11 The a-modification is a con- 
sistent choice for replacing the sine phase matching with a 
Gaussian phase matching, as it is implemented by modi- 
fying the width of the Gaussian and not by modifying the 
relation between K and the product ha ad hoc, as in (19). 



The second argument relays on a common criterium to 
choose a sensible value for a without recurring to ad hoc 
arguments. The common criterium is equal width at 1/e 
from the peak intensity [25^ . This criterium yields again a 
value of a 0.65. 



K 




0.1 0.2 0.5 1 2 



5 10 



Figure 10. K and i^az as functions of ha for (p = 0. 




Figure 11. K in green and the two different modifi- 
cations of the law ([5|. The dashed red line is the {a^ (3)- 
modification, the solid blue line is the a-only modification. 
The fits obviously match for 6cr < 1. 



(19) 
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Figure 12. Joint radial probability density of generating 
a pair of photons with 0AM = 10 and radial momentum 
qi and q2. This matrix is the matrix F for cp = used in 
the examples. 
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Figure 13. Joint radial probability density of generating 
a pair of photons with 0AM ^ = 10 and radial momentum 
Qi and q2. This matrix is the matrix F for = —4 used 
in the examples. 



4 Conclusions 



Appendix 

Here we report more in detail the function that we use to 
fill the matrix F. We start with the following two-photon 
wave function 



We developed a fast numerical approach to perform a 
numerical Schmidt decomposition without approximating 
the phase matching function, and for any phase mismatch. 
Thanks to the numerical approach we investigated the lim- 
itations of the Gaussian approximation and we proposed 
a scaling b ba, for a 0.65, that allows the Gaussian 
approximation with no phase mismatch to reproduce the 
results of the more complete sine phase matching. We also 
successfully verified the etendue argument by calculating 
K and Kaz for different values of the phase mismatch. 

The Schmidt modes that we found exhibit numerous 
features, the most remarkable being the grouping in local 
families of HG-like modes and the fact that they don't 
span more than one lobe on the diagonal of the matrices 
of joint radial probability density. We gave an interpreta- 
tion of these features by using the meaning of coherence 
length and the way of naturally assigning eigenmodes to 
potential wells. In fact, the other sets of HG-like radial 
modes are simply shifted by different amounts into other 
rings of the sinc-type SPDC emission. These other sets 
have similar HG-like radial structures but lower modal 
weights. The rather strict separation in discrete sets sug- 
gests that we might want to use the set number as third 
quantum number, i.e. to split the radial quantum number 
p in two quantum numbers instead. 



'^(qi,q2) = ATe-^l^^+^^l'sinc ( ^ + ^|qi -q2p 



(20) 



which can be written in cylindrical coordinates: 

/ 

= A/'exp 



■^{Ql + 12) - ^9l92 C0S((/?1 - ip2) 



\ 



X smc 



L 



(21(22 COs((^l - (/?2) 



V 




(21) 

= A/''exp(— acos((/9i — (/p2))sinc(/3 — jcos{(pi — (p2)) 

(22) 

where A/"' = A/'exp(-^f (9? + q^)) 

Similarly to what was done for the exp phase match- 
ing, we begin the evaluation of the Fourier transform with 
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ex.p{—ii{ipi — ^2)) as the Fourier kernel, which will give 
the eigenf unctions of a state of 0AM zb£. 

^^F{qi,q2) = 2TTN' / e-^"°^^sinc(/3 - 7COS(^)e-^^'^ 
^0 



Hi 



We have 



Hi= e-^^°^^^)sinc(/3 - 7cos((/p))e-^^^ dcp 
Jo 



(23) 



(24) 



1 ^it f 

^ J-L/2 Jo 

L 



27r 



I 

L/2 

( 

L/2 



'I\i\{-a 



{a+^t)cos(<p)-ii<f^^^^ (25) 

2^7 



-t) dt 



ijt) dt 



(26) 



(27) 



This last form of the integral can be approximated if 
the product ba can be considered sufficiently small. The 
amount of precision sought restricts the maximum value 
of Lji, but at Lji < 0.1 (which corresponds to ba < 0.2) 
the maximum error is less than 1%. The approximating 

function of ^^1+^2)^7^ is found to be 

f{qi,q2) = 27re-(^i+^2)/,(2gig2) sine [^/2LR {ql + ql) + ^ 

(28)" 

A matrix F is filled with entries {i^j) = Sy/Ijf{is,js), 
where 5 is the step, until the size N of the matrix is 
reached. Then = FF^ is used in the singular value de- 
composition. The resulting eigenvalues are Xi^p (at fixed 
i); the resulting eigenmode are (j)i^p{q) = y/qui^p^q). 

At sufficiently large x, the modifi ed Bessel function 
can be expanded as /^(x) ^ I\j2i\x. By combining this 
expansion with the mentioned multiplication by ^Jq\q2^ 
we find that the matrix F has a similar appearance as 
the original two-photon field of Eq. (4) , but with the roles 
of ((7i — ^2) and (gi + ^2) interchanged and a removal of 
the vector character. After Taylor expansion of the sinc- 
function, which works best in the considered limit Lji « 
1 for < 0, the resulting bi-exponential function yields 
a set of Hermit e- Gaussian radial Schmidt modes centered 
around a displaced maximum. 
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